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ABSTRACT 

V  , 

'An  analysis  of  seismic  propagation  in  the  crust  and  upper 
mantle  of  the  Kuriles/Kamchatka  region  has  been  completed  in 
support  of  the  MSS  site  evaluation.  Both  three-dimensional  ray 
tracing,  techniques  and  generalized  ray  theory  methods  have  beer, 
utilized  to  model  the  propagation  through  a  highly  complex, 
laterally  varying  structure.  The  structural  model  of  the  upper 
mantle  ar.d  descending  lithospheric  slab  developed  by  Veith 
(1974)  was  modified  for  this  study  to  include  the  observed 
change  in  angle  of  descent  of  the  slab  along  the  trench.  The 
intensity  of  incident ' seismic  energy  at  the  surface  was  mapped 
and  contoured  for  seismic  rays  refracting  through  the  three- 
dimensional  mantle  structure.  Synthetic  seismogram  profiles 
were  also  computed  to  investigate  the  waveguide  effects  of  the 
lithosphere  at  the  surface  and  to  analyze  the  contribution  of 
diffracted  and  multiple  reflected  energy,  particularly  within 
areas  of  classical  shadow  zones. 
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INTRODUCTION 

In  many,  if  not  most,  situations  involving  the  calculation 
of  seismic  propagation  within  the  earth,  the  assumption  that 
seismic  properties  vary  only  with  depth  may  be  used  successfully. 

Under  these  circumstances,  the  normal  vectors  to  the  propagating 
wavefronts,  or  rays,  are  well-behaved,  remaining  within  their 
initial  planes  of  propagation  through  the  earth.  As  a  result 
of  the  spherical  symmetry,  the  travel-times  and  amplitudes  of 
Seismic  phases  depend  only  on  distance  from  the  source  and  not  upon 
the  source-receiver  azimuth.  However,  in  the  presence  of  large- 
scale,  localized  heterogeneities  in  the  crust  or  upper  mantle, 
symmetry  assumptions  break  down  and  relatively  simple  techniques 
for  calculating  ray  paths,  travel-times,  and  amplitudes  are  inadequate 
This  situation  is  particularly  acute  when  the  lateral  heterogeneities 
are  in  the  near-source  region.  Complex  computational  schemes 
must  then  be  employed  to  follow  the  seismic  rays  point  by  point 
through  the  earth.  The  region  encompassing  northernmost  Japan, 
the  Kurile  Islands,  and  the  Kamchatka  peninsula  is  typified  by 
laterally  heterogeneous  crust  and  upper  mantle  structures.  Any 
detailed  evaluation  of  seismic  propagation  within  this  region 
must,  therefore,  take  into  account  these  variations  in  seismic 
properties . 

Planning  is  currently  underway  for  the  deployment  of  a 
down-hole,  ocean  bottom  seismometer  system,  the  Marine  Seismic 
System  (MSS) ,  roughly  several  hundred  kilometers  east  or  south- 
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east  of  the  Kurile  Islands.  Such  a  sophisticated  seismic 
monitoring  system  will  be  invaluable  in  increasing  not 
only  our  capabilities  for  nuclear  test  monitoring  and 
verification  but  also  our  knowledge  and  understanding  of 
the  earthquakes  and  tectonics  of  the  northwest  Pacific  and 
the  interaction  between  the  Pacific  and  Asian  plates  at  the 
Kurile  Trench.  Among  the  many  considerations  involved  in  the 
site  selection  for  the  MSS,  the  characteristics  of  seismic 
propagation  between  the  centers  of  seismic  activity  in  the 
.Kuriles/Kamchatka  region  and  the  proposed  MSS  sites  must  be 
considered.  This  report  examines  seismic  propagation  in  the 
Kuriles/Kamchatka  region  and  analyzes  the  potential  problems 
and  advantages  of  the  sites  under  consideration  from  a  seis¬ 
mologies!  viewpoint.  This  analysis  has  included  examination  of 
the  laterally  heterogeneous  structure,  distribution  of  seismicity 
types  and  spatial  distribution  of  source  mechanisms,  and  the 
propagation  characteristics  of  the  area,  including  both  geometric 
ray  propagation  and  non-geometric  paths  (diffracted  arrivals, 
head  waves,  etc.). 
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EARTH  MODEL  SELECTION 

t  The  Kuriles/Kamchatka  region  is,  geomorphologically ,  a 

typical  oceanic  trench/island  arc  structure.  To  the  east  of 
the  islands  is  a  deep  oceanic  trench.  The  islands,  as  well 
as.  much  of  the  Kamchatka  peninsula,  are  highly  volcanic. 
Seismically,  the  area  is  extremely  active.  Hypocentral 
depths  range  from  the  surface  to  600-700  km  in  the  mantle. 

This  kind  of  geomorphologic  and  tectonic  feature  is  very 
successfully  described  within  the  framework  of  seafloor 
spreading  and  plate  tectonics. 

The  theory  of  plate  tectonics  is  now  widely  accepted  by 
earth  scientists  and  is  regularly  employed  in  interpreting 
geological- and  geophysical  observations  and  explaining  the 
tectonic  environment  and  evolution  of  the  earth.  The  prime 
concept  of  plate  tectonics  is  that  the  outer  layer  of  the 
earth,  the  lithosphere,  is  broken  into  a  number  of  relatively 
rigid,  individual  plates,  each  moving  relative  to  one  another. 
The  driving  mechanism  for  this  motion  is  generally  presumed  to 
be  thermally-derived.  The  plates  interact  in  primarily  three 
different  modes.  In  regions  where  the  motion  is  divergent, 
new  material  upwells  from  deeper  in  the  mantle, forming  long 
ridge-like  or  rifting  structures.  The  mid-oceanic  ridge  system 
is  the  most  prominent  such  feature.  Earthquake  mechanisms 
associated  with  this  tectonic  activity  are  generally  tensional 
or  normal  faulting.  The  second  interaction  occurs  where  two 
plates  are  moving  in  roughly  parallel  directions,  a  transform 
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fault  type  of  relationship.  Such  regions  (e.g.,  the  San 
Andreas  fault  zone  in  California)  are  typified  by  strike-slip 
earthquake  faulting.  The  third  region  of  interaction  is  a 
convergent  zone,  the  type  of  zone  of  most  interest  to  us  in 
this  study,  in  a  cypical  convergent  zone,  one  plate  is  under¬ 
thrusting  the  second.  In  Figure  1,  from  Ringwood  1?69)  ,  the 
general  framework  of  this  zone  is  illustrated.  The  seismicity 
in  convergent  zones  ranges  in  depth  from  very  shallow  to  quite 
deep  (600-700  km)  following  a  dipping  zone  of  seismic  activity, 
called  a  Benioff  zone.  The  faulting  types  within  this  zone 
include  normal,  thrust,  and  strike-slip  mechanisms.  However, 
the  different  kinds  of  faulting  are  not  randomly  mixed  but 
occur  in  systematic  patterns  based  on  the  local  tectonic 
environment  (for  example,  in  the  region  of  maximum  bending  on 
the  top  of  the  underthrust  slab,  the  earthquakes  show  tensions! 
faulting).  xn  some  regions,  such  as  the  Northeast  Japan  arc, 
the  seismicity  falls  into  a  double  plane,  clearly  resolving 
the  downgoing  slab  (Figure  2) . 

The  Kurile  Island/Kamchatka  region  is  a  zone  of  con¬ 
vergence  between  the  Pacific  and  Asian  plates.  The  Pacific 
plate  is  sharply  underthrusting  the  Asian  plate  in  this  region. 
The  angle  of  underthrusting,  and  therefore,  the  dip  of  the 
Benioff  zone  defining  the  focus  of  seismic  activity,  is  steeply 
dipping.  It  is  this  highly  complex,  laterally  heterogeneous 
region  that  must  be  modeled  in  order  to  analyze  seismic  prop¬ 
agation  in  support  of  the  MSS  site  evaluation.  Fortunately, 
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Figure  1.  Generalized  framework  of  the  spreading  and  descending 
plates  (from  Ringwood,  1969). 
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Figure  2.  Two-layered  deep  seismic  plane  under  the  Northeast 
Japan  Arc  (from  Hasegawa,  Umino,  and  Takagi,  1976). 
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this  region  has  received  a  significant  amount  of  attention  in 
the  past  several  years,  and  an  appropriate  crust  and  upper 
mantle  model  can  be  constructed  with  some  degree  of  confidence. 

The  most  thorough  study  of  the  seismicity  and  structure 
of  the  Kurile  Islands  region  was  performed  by  Veith  (1974) . 

Veith  relocated  over  three  thousand  earthquakes  in  the  Kurile 
Islands  using  source-region/station  time  corrections.  This 
permitted  the  determination  of  very  precise  hypocentral  locations. 
Moreover,  useful  source  mechanisms  were  determined  for  over 
four  hundred  of  these  events.  This  study  was  thus  able  to 
very  clearly  delineate  not  only  both  the  non-laterally  heter¬ 
ogeneous  upper  mantle  structure  in  the  area  and  the  structure 
and  velocity  contrast  to  the  downgoing  slab  but  also  the 
spatial  variations  and  trends  in  faulting  mechanism. 

In  Figure  3,  Veith' s  (1974)  variation  of  seismic  velocity 
with  depth  in  both  the  slab  and  the  surrounding  mantle  is 
presented.  The  velocity  contrast  between  the  upper  mantle  and 
the  descending  plate  ranges  from  about  4%  at  600-700  km  depth 
to  a  maximum  of  approximately  14%  at  about  200  km  depth.  Such 
a  dramatic  contrast  does,  as  will  be  discussed  in  a  later 
section,  cause  fair.'.y  dramatic  alterations  of  ray  paths.  For 
shallow  (less  than  -100  km)  depth  events  in  the  Kuriles,  how¬ 
ever,  the  presence  of  slab  at  depth  is  of  little  or  no  con¬ 
sequence  for  propagation  to  the  proposed  MSS  sites.  The  most 
important  model  features  for  those  cases  is  the  rather  typical 
oceanic  structure  with  a  well-developed  lid  and  underlying  low 
velocity  zone.  This  lid  structure  provides  a  strong  wave-guide 


for  trapped  energy.  Those  effects  are  treated  separately 
later  in  this  report. 

The  slab  model  shown  in  Figure  3  is  Veith's  (1974)  com¬ 
posite  or  average  model,  the  dip  angle  of  the  descending  slab 
being  about  50°.  The  angle  of  descent  actually  varies  along 
the  trench  from  about  40°  north  cf  Hokkaido  to  about  60°  under 
Kamchatka  (Veith,  personal  communication).,  In  this  study,  the 
three-dimensional  slab  model  used  has  this  warped  configuration. 
In  addition,  the  length  of  the  slab,  and,  correspondingly  the 
maximum  depth  attained,  as  defined  by  the  observed  seismicity, 
is  also  a  function  of  position  along  the  trench.  In  the  shallow 
dipping  regions  at  the  southern  extreme,  the  maximum  depth  is 
350-400  kilometers;  in  the  north,  the  maximum  depth  is  about 
700  kilometers.  The  velocity  model  for  the  warped  slab,  how¬ 
ever,  still  corresponds  to  the  velocity  vs.  depth  profile  shown 
in  Figure  3. 
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APPROACH 

The  basic  approach  in  our  evaluation  of  the  potential  MSS 
sites  in  the  Kuriles  region  has  been  to  select  a  variety  of 
representative  hypocentral  locations  within  the  area  and  to 
then  calculate  an  appropriate  suite  of  the  propagation  paths 
from  those  hypccenters  to  the  region  encompassing  all  of  the 
possible  MSS  locations  east  of  the  trench  axis.  This  approach 
to  ray  propagation  in  the  region  was  the  best  means  to  system¬ 
atically  analyze  propagation  from  a  large  number  of  potential 
source  sites  to  a  fairly  broad  region  in  which  the  MSS  might 
be  located. 

Seismic  propagation  in  the  region  was  divided  into  two 
distinct  categories  and  different  techniques  were  employed 
for  each  class  of  problems.  The  first  grouping  contains  all 
cases  in  which  the  hypocentra’.  location  under  consideration 
was  deeper  than  100  km.  In  such  cases,  the  presence  of  the 
descending  slab  exerted  a  dominent  propagation  effect.  Seismic 
propagation  was  modeled  using  solely  three-dimensional  ray¬ 
tracing  methods  since  such  deep  sources  are  extremely  inefficient 
in  pumping  trapped  energy  into  the  shallow  lithospheric  wave¬ 
guides.  The  second  category  of  problems  involved  locations  at 
depths  less  than  100  km  especially  very  shallow  sources.  These 
represent  a  very  different  class  of  propagational  problems. 

At  shallow  source  depths,  a  significant  proportion  of  the  seismic 
energy  is  trapped  within  the  relatively  efficient  lithospheric 
waveguides.  The  behavior  of  such  radiation  is  independent  of 
the  deep,  downgoing  slab  effects,  but  is  very  dependent  upon 
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frequency,  fine  structure  (especially  sediment  thicknesses, 
and  source  depth).  Multiply-reflected  rays,  head  waves,  and 
other  non-geometric  paths  proliferate.  We  have  employed 
complete  generalized  ray  techniques  to  model  this  kind  of 
propagation.  Because  of  the  sensitivity  to  fine  structure, 
however,  this  part  of  the  analysis  has  been  conducted  in  only 
generic  form.  More  detailed  analysis  can  be  employed  for  a 
specific  path  or  suite  of  paths  once  a  site  has  been  selected. 
Three-dimensional  ray-tracing  analyses  have  also  been  conducted 
and  some  examples  are  presented  for  shallow  sources.  However,  as 
will  be  seen,  the  basic  oceanic  structure  in  the  region  causes 
direct  rays  to  be  largely  channeled  out  of  the  region  of  interest. 

The  propagational  modeling  presented  here  also  includes 
the  effects  of  differing  radiation  patterns.  The  study  by 
Veith  (1974)  of  this  region  was  used  to  determine  the  dominant 
faulting  types  throughout  the  area  and  the  resulting  radiation 
pattern  effects  were  convolved  into  our  analysis.  In  Veith' s 
study,  five  distinct  types  of  faulting  were  found;  transverse, 
shallow  normal,  thrust,  compressional,  and  tensional.  Trans¬ 
verse  mechanisms  occur  at  all  depths  within  the  plate  and  allow 
for  both  the  differential  rate  of  underthrusting  along  the 
trench  (7.5  cm/yr  in  the  north  to  8.5  cm/yr  in  the  south; 

Le  Pichon,  1968)  and  for  the  warping  or  uneven  bending  of  the 
slab.  Such  mechanisms  comprise  about  20%  of  the  faulting  over¬ 
all  and  greater  than  30%  of  events  deeper  than  70  kilometers. 

The  events  are  nearly  dip-slip  slip  with  the  northern  segment 
moving  down  relative  to  the  southern  side.  Shallow  normal 
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faulting  generally  occurs  only  in  areas  of  the  oceanic  trench. 

The  third  type  of  faulting,  thrust  faulting,  is  the  dominent 
mode  observed  in  the  region  comprising  about  45%  of  all  observed 
mechanisms.  The  final  two  types  of  earthquakes  occur  in  rel¬ 
atively  narrow  zones  at  depth  and  it  has  been  suggested  that 
they  may  be  related  to  a  phase  transition  in  the  downgoing 
slab  (Veith,  1974) . ' 

Figure  4  shows  the  five  focal  mechanisms  which  Veith  found 
characteristic  of  the  Kurile  region.  Shaded  and  unshaded  areas 
of  the  hemispheres  correspond  to  the  polarity  of  tne  propagating 
wavefront.  Rays  leaving  the  source  through  the  shaded  regions 
have  compressional  first  motions  while  those,  leaving  through 
unshaded  areas  have  dilatational  wavefronts.  The  energy  radiated 
from  an  earthquake  source  is  distributed  according  to  the  rad¬ 
iation  pattern.  No  energy  escapes  the  source  along  nodal  lines 
in  the  focal  mechanism  while  maximum  energy  is  radiated  from 
the  center  of  each  region.  Notice  the  similarity  in  fault 
orientation  between  the  normal  and  thrust  mechanisms  and  between 
compressional  and  tensional.  These  pairs  essentially  differ 
only  by  wavefront  polarity  which  is  not  a  necessary  parameter 
in  the  calculation  of  radiated  energy.  Fcr  this  reason  we 
only  consider  three  basic  fault  mechanisms  in  the  Kuriles, 
with  the  understanding  that  first  motions  may  be  either  compress- 
iotial  or  dilatational  for  two  of  these  mechanisms.  . 

The  most  useful  end  product  of  the  ray  tracing  modeling 
are  spatial  contours  delineating  the  intensity  of  incident 


direct  arrivals,  referred  to  as  incident  energy-density  contours. 
The  generalized  ray  calculations  yield  synthetic  seismogram 
contours.  Both  provide  useful  tools  for  evaluating  the  suit¬ 
ability  of  potential  station  site'. 

RAY  TRACING 

Three-dimensional  seismic  ray  tracing  techniques  have 
been  applied  in  several  studies  of  island  arc  regions.  Jacob 
(1972)  and  Davies  and  Julian  (1972)  sought  to  explain  anomalous 
travel  time  residuals  from  the  underground  nuclear  explosion 
Longshot  fired  in  the  Aleutians.  Large  lateral  variations  in 
the  upper  mantle  caused  by  a  descending  slab  of  lithospheric 
material  severely  distorted  raypaths  from  this  event.  In 
these  studies,  as  in  oiirs,  a  regional  earth  model  was  adapted 
for  use  in  a  three-dimensional  ray  tracing  code.  Rays  were 
traced  through  the  slab  model  and  travel-time  residuals  compared 
to  those  observed.  In  another  study  (ToksSz,  et  al.,  1971), 
ray  travel  times  for  the  Tonga-Fiji  region  were  correlated 
with  three-dimensional  ray  tracing  results  and  theoretical 
travel  times  based  on  calculated  temperature  fields.  More 
recent  work  by  Engdahl  (1980),  has  produced  a  highly  detailed 
earth  model  of  the  Aleutian  island  arc  region.  Large  travel¬ 
time  data  sets  were  used  in  conjunction  with  a  three-dimensional 
ray  tracing  code  to  accurately  determine  the  model. 

Ray  theory  is  particularly  applicable  to  wave  propagation 
problems  involving  inhomogeneous  media  with  large  linear 
dimensions  compared  to  the  dominant  wavelengths.  In  such  cases. 
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motion  near  a  wavefront  in  examined  and  propagation  r.ay  be 
characterized  implicitly  by  an  expanding  surface  T(X,E),  where 
T  is  the  travel  time  between  the  source  at  E  and  a  point  X 
on  the  wavefront.  The  possible  wavefronts  are  constrained  by 
requiring  that  T  satisfy  the  eikonal  equation: 

(VT)2  =  1/C2  (1) 

where  C  equals  the  local  P  or  S  wave  velocity  along  the  wave- 
front.  Since  we  are  interested  in  tracing  the  surface  explicitly 
as  a  function  of  travel  time,  the  form  X  =  X(T)  is  more  useful. 
Ray  vectors,  S,  are  defined  as  normals  to  the  wavefront,  equiv¬ 
alent  to  VT  with  magnitudes  equal  to  inverse  local  velocity 
(1/C (X) )  or  N  (X) ,  the  ray  slowness.  From  these  definitions, 
a  set  of  differential  equations  describing  ray  propagatio.  may 
be  expressed  in  Cartesian  coordinates  as  follows : 


dX./dT  =  S./N2 

l  i 

(2) 

dlT/dT  =  (1/NJ  dN/d^ 

(3) 

The  ray  tracing  code  solves  these  equations  for  X^  and  in 
spherical  coordinates  by  standard  numerical  methods. 

The  minimum  requirement  for  a  solution  to  the  ray  tracing 
equations  is  that  first  derivatives  of  N(X)  be  continuous  at 
each  point  X.  A  cubic  spline  interpolation  scheme  was  adapted 
for  this  purpose  since  by  definition  all  first  derivatives  at 
X  are  forced  to  be  continuous.  Using  this  scheme,  cubic  poly- 


nomials  are  fit  between  sampled  N(X)  values  as  ray  propagation 
progresses  and  spatial  derivatives  are  calculated.  This  pro¬ 
cedure  poses  a  problem  however,  for  accurately  representing  a 
layered  earth  model. 

Polynomial  approximations  such  as  least  squares  are 
notorious  for  their  resistance  to  bending- around  sharp  corners, 
but  spline  approximation  is  especially , unyielding.  Step  dis¬ 
continuities  in  the  data  are  splined  together  with  continuous 
cubic  polynomials  producing  numerous  overshoots  and  undershoots. 
Layered  earth  models  are  especially  prone  to  this  error,  being 
discrete  functions  with  numerous  discontinuous  velocity  steps. 
Figure  5  illustrates  this  point.  Veith's  original  upper  mantle 
earth  model  (homogeneous  part)  is  shown  here  with  a  spline 
approximation  curve.  The  large  irregularities  in  the  spline 
would  noticeably  alter  ray  propagation  paths  based  on  Veith’s 
model. 

Since  the  ray  tracing  code  uses  spline  interpolation, 
the  initial  model  shown  in  Figure  5  must  be  altered  in  such  a 
way  that  when  it  is  splined,  the  result  matches  the  original 
model.  Initially,  the  earth  model  is  sampled  evenly  in 
radius  (depth)  every  25.93  km  to  a  depth  of  1200  km.  This 
is  a  necessary  formating  step  to  prepare  the  model  for  input. 
Note  that  beyond  the  lower  depth  limit  of  Veith's  upper  mantle 
model  an  Anderson  and  Hart  (1975)  earth  model  is  used.  Ideally 
of  course,  the  spline  problem  would  vanish  as  the  sampling 
increment  decreased.  The  computational  cost  involved  was  tco 
high  to  justify  a  densely  sampled  earth  model.  The  result  of 
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discrete  sampling  is  shown  in  Figure  6a  and  is  almost  identical 
to  the  original  as  one  would  expect.  Next,  in  Figure  6b,  the 
discretely  sampled  velocity  structure  is  varied  to  prevent 
overshoots  and  undershoots  present  in  Figure  6a.  This  spline 
is  now  representative  of  the  original  model  as  shown  in  Figure  6c 
The  modified  velocity  profile  in  6b  is  suitable  for  input  to  the 
code  since  its  spline  interpolation  closely  approximates  the 
initial  model. 

So  far  we  have  only,  addressed  the  problem  of  preparing  the 
nomogeneous  earth  model.  Heterogeneous  structure  is  input  as 
a  perturbation  to  the  layered  earch  model  and  a  similar  pro¬ 
cedure  is  used  to  integrate  the  perturbation  into  the  overall 
model.  In  this  case  the  perturbations  are  modified  slightly 
such  that  when  added  to  the  standard  earth  model  the  final 
profile  is  well-represented  by  a  spline  interpolation.  Figure  7 
shows  the  sampling  grid  used  to  input  velocity  perturbations 
due  to  the  descending  lithospheric  slab.  This  portion  of  the 
model  is  sampled  evenly  in  R  and  9,  where  dR  is  the  same  as 
before  and  d8  *  0.1214  degrees.  Perturbations  are  tapered  to 
zero  exponentially  outside  the  slab  and  attain  a  maximum  of 
1.2  km/sec  near  35C  km  depth.  The  tapering  was  introduced  to 
help  reduce  splining  errors. 

The  ray  tracing  i3  performed  in  spherical  coordinates  (R, 

9,  <> )  with  respect  to  a  specified  auxiliary  coordinate  system 
to  aid  in  the  placement  of  complicated  structures.  Inhomogeneous 
media  are  input  in  R,  9  coordinates  as  perturbations  of  a 
homogeneous  earth  model  C  (X) .  These  may  be  discrete  or 
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Figure  6.  The  evenly  sampled  (dR£-25Km)  version  of  Veith's  velocity  profile  is  shown  in  6a 
(dotted  line)  next  to  the  original.  Fig.  6b  illustrates  how  small  modifications 
to  this  model  prevent  overshoots  from  splinning.  Finally,  the  spline  curve  from 

6b  is  shown  in  6c  overlaying  the  original  model.  Since  the  ray  tracing  code  use 
saline  interDolation .  the  modified  profile  (6b)  is  incut. 


KURILE  VELOCITY  MODEL 
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functional  representations.  The  perturbed  two-dimensional 
region  is  then  extended  along  the  $  axis  and  tapered  at 
specified  end  points  to  generate  a  three-dimensional  structure. 
The  region  must  therefore  have  a  standard  cross-section  with 
respect  to  the  $  coordinate. 

A  generalized  coordinated  system  is  initially  defined  by 
Eulerian  rotation  angles  (e.g. .,  Goldstein,  1950) .  The  Euler 
angles  0  ,  <p  ,  and  ip  specify  a  series  of  three  rotations  from 

r  r  r 

an  initial  system  of  axes,  X,  Y,  Z,  to  a  new  system,  X  ,  Y  ,  Z  , 

P  P  P 

and  should  not  be  confused  with  the  spherical  coordinates  R,  9, 

and  <j>.  In  this  case  the  initial  axes  correspond  to  geographic 

coordinates  with  the  positive  Z  axis  at  90°  North  latitude 

(north  pole)  and  the  positive  X  axis  at  0°  longitude  on  the 

equator.  e  is  then  the  angular  distance  between  the  Z  and  Z 
P  P 

axes.  4>p  specifies  a  counterclockwise  rotation  of  the  X  axis 

to  a  i.ew  position  E  in  the  X-Y  plane,  ip  is  the  angular  distance 

between  X  and  E  subtended  from  a  counterclockwise  rotation  in 
P 

the  Xp“Yp  Plane.  Figure  8  shows  the  generalized  coordinate 
system  used  to  model  the  Kurile  subdue tion  zone. 

Veith's  velocity  profile  of  the  Kurile  region  is  represented 
as  described  previously,  in  discrete  form  as  a  function  of  R 
and  0.  Since  this  is  an  average  model  it  is  valid  as  a  first 
approximation  at  every  point  along  the  arc.  Choosing  the 
trench  axis  cs  a  reference  line,  we  would  like  the  coordinate 
9  to  be  constant  along  this  line  in  order  for  any  cross-section 
of  the  region  to  have  identical  velocity  profiles  when  the  moael 
is  extended  to  three  dimensions.  This  places  a  constraint  cn 
the  auxiliary  pole  location  Zp. 


Figure  8. 


Ray  tracing  is  performed  in  spherical  coordinates 
(R,  9,  <p)  with  respect  to  an  auxiliary  coordinate 
system  (X  ,  Y  ,  Z  ) .  Rotation  from  a  geographic 
system  (X?  Y,^Z)  ?s  specified  by  the  Euler  angles 
0  ,  <J>  ,  and  .  The  Kurile-Kamchatka  arc  is 
p?ace8  into  tftis  auxiliary  system  by  requiring 
the  pole  Z  to  be  equidistant  from  all  points  on 
the  arc.  ?n  addition,  the  prime  meridian  (4>=0  line) 
is  oriented  to  intersect  the  arc's  southern-most 
extension. 


The  Kurile-Kamchatka  arc  was  placed  into  the  auxiliary 
coordinate  system  by  requiring  the  pole  Zp  to  be  equidistant 
from  all  points  on  the  trench  axis.  To  accomplish  this,  the 
angular  distance  to  the  pole  was  varied  until  a  best  fit  to 
the  trench  axis  was  found.  Davies  and  Julian  (1972)  placed 
the  Aleutian  arc  by  the  same  technique  finding  that  10 6  was 
an  optimal  polar  angle.  In  our  case,  a  polar  distance  of  35° 
fit  the  data  to  within  an  average  error  of  0.10°.  The  following 
is  a  list  of  7  trench  axis  points  and  their  corresponding  angle 
to  the  auxiliary  pole  Zp: 


Latitude 

Longitude 

Angle 

42. 7N 

148. 4E 

34.98 

43. 6N 

150. 0E 

35.07 

•  45. ON 

152. 4E 

35.18 

47. N 

155. 0E 

35.00 

48. 6N 

157. 0E 

34.84 

50. 5N 

160. 0E 

34.98 

52. 8N 

162. 4E 

34.60 

(reference) 


The  actual  location  of  the  Z^  axis  (pole)  in  geographic  coordinates 
is  61.18°  N  and  97.86°E.  The  prime  meridian  in  this  coordinate 


system  intersects  the  southern  most  point  on  the  Kurile  trench 


axis.  Hence,  the  azimuthal  angle  $  equals  0.°  at  this  point. 

Since  the  Kurile  arc  is  approximately  25°  long  and  8.5°  wide, 
the  regional  model  is  defined  for  0  =  26.5°  to  35°  and  <J>  =  0°-24°. 
The  maximum  depth  of  the  lateral  heterogeneity  is  700  km. 

The  velocity  profile  determined  by  Veith  (Figure  7)  is  an 
average  model  of  the  region.  In  reality,  the  dip  angle  of  the 
slab  is  known  to  vary  by  20°  from  north  to  south.  In  order  to 
incorporate  this  feature  into  our  model,  the  velocity  profile 
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of  the  slab  was  made  functionally  dependent  or  R  and  <j>.  At  the 
extreme  southern  end  of  the  trench,  where  <p  -  0.°,  the  dip  angle 
was  set  to  40°  by  shifting  the  original  perturbations.  The 
shifting  function  varies  dip  angle  continuously  up  to  60°  at 
the  northern  most  point  of  the  trench.  Northern  and  southern 
cross-sections  of  the  warped  slab  model  are  shown  in  Figure  9. 
Velocity  contour  lines  are  drawn  every  0.5  km/sec  from  11.50  km/sec 
just  above  1200  km.  The  obvious  area  of  discontinuity  outlined 
here  corresponds  to  the  descending  slab. 

Full  three-dimensional  variance  of  the  model  was  accomplished 
by  functional  manipulations  of  a  standard  cross-section.  This 
technique  may  be  applied  to  almost  any  complicated  structure 
given  some  fundamental  profile.  The  ability  to  rotate  the 
coordinate  system  also  aids  in  placing  structures  as  we  have 
shown  for  the  Kurile  arc.  Together,  these  two  features  make 
this  ray  tracing  code  a  powerful  tool  in  geophysical  research. 

To  fully  appreciate  the  three-dimensional  nature  of  the 
slab  models,  a  radial  velocity  contour  map  was  constructed. 

Figure  10  shows  the  regular  and  warped  slab  models  respectively 
at  four  depths.  Each  contour  separates  differing  velocity  regions, 
thus  outlining  the  descending  slab. 

Utilizing  the  ray  tracing  techniques  described  above, 
incident  energy  density  contour  maps  have  been  constructed  in 
the  Kurile  region  based  on  data  collected  from  rays  propagated 
through  the  warped  slab  model.  Sources  for  these  rays  were 
selected  to  sample  local  seismicity  data  as  previously  discussed 
and  the  three  basic  focal  mechanisms  selected  were  thrust. 


Figure  10.  Radial  cross-sections  of  the  regular  (dotted  line)  and  warped  slab  models 
at  four  depths.  Velocity  contours  outline  the  slab's  position,  varying 
from  8.0  km/sec  at  115  km  to  10.5  km/sec  at  640  km. 
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transverse,  and  compressional.  Relative  energy  of  a  seismic 
ray  depends  significantly  on  both  source  radiation  pattern  and 
geometrical  spreading. 

In  a  homogeneous  earth,  the  geometrical  spreading  correction 
is  a  simple  function  G(A)  where  A  is  the  angular  distance  between 
the  epicenter  and  ray  emergence  point.  The  structural  complexity 
of  the  Kurile  subduction  zone  however,  produces  anomalous  focusing 
and  defocusing  of  rays  which  distort  normal  spreading.  Corres¬ 
ponding  irregularities  in  energy  density  (E(A))  are  produced 
at  the  surface  as  a  function  of  both  A  and  azimuth.  In  fact, 

E (A)  varies  proportionally  to  f_^/6A  where  e^  is  the  emergence 
angle  at  the  source.  The  exact  relation  is  given  by  Bullen  (1965): 

E  (A)  =  (I/Rq2)  *  (cos  (e.^  )  /sin  (A)  *sin  (eQ) )  *  |  6ei/(SA  |  (4) 

where  Rq  equals  the  earth's  radius,  I  is  the  energy  per  unit 
solid  angle  emitted  from  the  source,  eQ  is  the  emergence  angle 
at  the  surface,  and  A  is  the  angular  distance  between  the  epi¬ 
center  and  ray  emergence  point.  E(A)  is  then  the  energy  per 
unit  area  of  the  portion  of  the  wavefront  emerging  at  angle  eQ. 

For  computational  ease,  the  derivative  5e^/6A  is  determined 
discretely. 

In  our  modeling,  rays  were  emitted  from  the  sources  at 
a  constant  increment  of  take-off  angle  (usually  2.00°),  hence 
<$e^  is  a  constant.  As  a  result,  the  energy  relation  is  most 
dependent  on  6 A  which  is  always  changing.  Where  ray  focusing 
occurs  the  change  in  A  is  small  thus  raising  the  incident  energy- 
density.  These  areas  are  known  as  "bright  spots"  or  caustics 
in  the  case  of  perfect  focusing. 
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Having  once  specified  a  source  mechanism  by  the  fault  plane 
orientation  (dip,  strike,  and  rake),  a  quadrupole  radiation 
pattern  is  determined.  Rays  exiting  the  source  have  relative 
energies  corresponding  to  their  angle  and  azimuth  of  departure. 
Energy  variation  from  radiation  pattern  is  a  straightforward 
calculation  (Langston  and  Helmberger,  1975)  and  the  resulting  . 
value  for  each  ray  may  simply  be  multiplied  by  E(A).  Radiation 
values  vary  between  a  minimum  of  0,  for  rays  exiting  along 
nodal  planes,  and  a  maximum  of  inverse  source  velocity  squared. 
Polarity  of  the  wavefront  is  indicated  by  the  sign  of  the 
radiation  value. 

Extensions  to  the  ray  tracing  code  originally  developed 
by  Bruce  Julian  in  1974  make  it  possible  to  calculate  regional 
energy  density  given  an  appropriate  earth  model.  In  this  case 
Veith's  model  of  the  Kuriles  is  adapted  for  input  to  the  code 
and  surface  energy  density  maps  constructed  to  the  southeast 
of  the  island  arc. 

Five  source  regions,  coded  KAn-KEn,  were  examined  at  varying 
depths  along  the  slab.  Depths  are  numerically  coded  with  n  =  1 
corresponding  to  the  shallowest  source  at  150  km,  and  n  =  2,  3, 
and  4  corresponding  to  progressively  deeper  hypocentral  locations 
Figure  11  illustrates  the  epicentral  locations  of  the  source 
regions  used  in  this  phase  as  well  as  the  four  proposed  MSS 
sites.  Rays  from  these  sources  were  traced  into  a  region  which 
enclosed  the  four  proposed  OBS  sites  by  specifying  hypocenter 
location  and  propagation  direction.  More  specifically,  rays 
were  shot  from  each  of  16  total  sources  in  2.0®  degree  intervals 
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of  take-off  angle  and  4.0°  intervals  of  azimuth  into  this  area. 
In  most  cases  at  least  400  rays  were  necessary  for  adequate 
coverage  of  the  region. 

Since  the  earth  model  is  fully  three-dimensional,  rays 
were  not  just  refracted  within  their  initial  departure  plane. 
Many  raypaths  were  forced  out  of  the  plane  by  the  velocity- 
structure  producing  unique  focusing  effects  at  the  surface, 
figures  12a  and  12b  are  cross-sections  of  the  earth  model 
showing  numerous  raypaths  from  sources  KD4  and  KC2  respectively. 
Rays  bending  out  of  the  plane  here  are  projected  into  the 
initial  ray  departure  plane.  Tick  marks  on  the  surface  are  at 
1°  intervals  and  the  depth  is  given  in  kilometers.  An  outline 
of  the  slab  is  produced  by  plotting  velocity  contours  every 
0.5  km/sec  as  was  done  in  Figure  9.  Note  the  increasing 
distortion  in  the  raypath  as  the  rays  strike  the  slab-earth 
interface  at  increasingly  oblique  angles  in  12a.  The  refraction 
of  these  rays  is  governed  pimply  by  Snell's  Law.  In  the  case 
of  12b  we  see  that  a  low  velocity  zone  at  150  km  is  actually 
trapping  rays  by  bending  them  away  from  the  interface  as  they 
strike  a  higher  velocity  medium  above.  Raypaths  in  Figure  12b 
are  severely  distorted  by  the  structure  causing  shadow  zone 
formation  just  beyond  the  trench  axis. 

As  source  locations  approach  the  surface,  the  low  velocity 
zone  traps  an  increasing  number  of  rays,  producing  an  8°  wide 
shadow  zone  southeast  of  the  trench  axis.  Figure  13  is  a 
cross  section  of  these  raypaths  from  a  surface  focus  event. 

The  take-off  angle  increments  have  been  made  very  small  in 
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Cross-sectional  map  of  raypaths  from  a  surface 
focus  event  in  Kamchatka.  A  10°  shadow  is  produced 
by  an  upper  mantle  low  velocity  zone. 
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order  to  show  more  detail.  The  implication,  of  course,  is 
that  direct  seismic  arrivals  from  very  shallow  sources  will  be 
low  amplitude,  long  period  signals  at  any  of  the  proposed  sites. 

A  treatment  of  this  problem  using  full  Cagniard  methods  is 
presented  later. 

Figures  14a  and  14b  present  one  mode  of  analyzing  the 
ray  tracing  results.  In  these  figures,  the  locus  of  emergence 
of  each  ray  is  plotted  as  a  small  vector.  The  vectors  indicate 
the  apparent  back  azimuth  of  each  ray  as  determined  by  the  angle 
of  emergence.  There  are  several  interesting  features  of  the 
propagation  which  become  apparent  in  these  figures.  In  Figure  14a, 
the  hypocenter  is  located  at  a  depth  of  600  km  near  the  northern 
terminus  of  the  trench.  The  radical  bending  of  some  rays 
from  this  location  can  be  seen  rather  dramatically  in  the  area 
directly  east  of  the  epicenter  where  the  apparent  back  azimuth 
is  rotated  as  much  as  180°  from  the  true  back  azimuth.  In 
Figure  14b,  another  effect,  that  of  horizontal  caustics,  is 
illustrated  from  an  event  near  the  southern  end  of  the  trench. 

This  caustic  can  be  seen  in  the  curving,  increasingly  dense 
pattern  of  arrivals  from  the  southeast. 

Although  the  presentation  in  Figure  14  can  yield  interesting 
and  important  information,  it  suffers  the  limitation  of  being 
a  discrete,  pointwise  representation  of  the  data.  This  limit¬ 
ation  has  been  surmounted  through  the  implementation  of  the 
incident  energy-density  contouring  described  earlier.  Figures  15 
through  30  show  the  incident  energy-density  distributions  from 
the  16  hypocentral  locations  shown  in  Figure  12.  Each  figure 


Figure  14.  Surface  maps  of  ray  emergence  points  for  two  Kurile  sources.  Arrows 
indicate  emergence  azimuth.  Structural  interaction  anomalously  bends 
rays  around  and  produces  unusual  focusing  effects. 
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in  this  group  is  a  composite  of  four  individual  plots.  In 
the  lower  left,  a  back-azimuth  vector  plot  of  the  type  shown 
in  Figure  14  is  displayed  for  ah  explosion  source  (no  radiation 
pattern).  In  the  upper  left,  upper  right,  and  lower  right 
quadrants  of  each  figure  are  the  energy-density  plots  for  the 
same  hypocentral  location  but  with  the  radiation  patterns  con¬ 
volved  into  the  result.  These  radiation  patterns  are  indicated 
by  the  focal  spheres  illustrated  in  each  plot.  For  clarity, 
the  contours  are  illustrated  in  color.  Changes  in  color  from 
blue  to  green  to  red  to  yellow  indicate  decade  changes  in 
incident  energy.  Blank  or  uncolored  zones  indicate  either 
complete  shadow  zones  or  regions  not  included  in  the  calculation 
The  proposed  MSS  sites  as  well  as  the  locus  of  the  Kamchatka 
peninsula  and  Kurile  Islands  have  been  plotted  on  a  transparent 
overlay  for  use  with  these  figures  (see  pocket,  inside  back 
cover).  In  all  cases,  the  propagational  effects  due  to  the 
downgoing  slab  dominate  the  broad  features  of  the  incident 
energy.  The  detailed  distribution  of  relative  energy  is, 
however,  more  sensitive  to  the  assumed  radiation  pattern  of 
the  source.  It  should  be  noted  that  the  energy  density  levels 
are  not  normalized  by  figure  but  rather  all  figures  are  plotted 
on  the  same  intensity  scale. 

It  is  also  possible  to  utilize  this  ray  tracing  approach 
to  analyze  potential  problems  m  source  mechanism  studies. 

Due  to  the  extreme  lateral  refractions  in  this  structure, 
nodal  planes  are  often  not  clearly  resolved.  This  result  is 
implied  from  Figures  15-30.  In  these  figures,  it  is  at  best 
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rare  to  observe  absolute  null  regions  corresponding  to  radiation 
pattern  nulls.  In  Figure  31,  this  result  is  illustrated  more 
directly.  Figures  31a  and  31b  plot  the  emergence  points  and 
radiation  polarity  from  two  deep  thrust  events  at  the  northern 
end  of  the  trench,  the  first  at  a  depth  of  325  km,  the  second 
at  550  km.  Compressional  and  dilatational  arrivals  are 
represented  by  open  circles  and  crosses  respectively.  The 
strong  overlap  of  the  two  polarities  clearly  illustrate  the 
problem.  This  result  is  presented  not  as  a  necessarily 
important  siting  consideration  but  rather  as  an  indication 
of  potential  ambiguities  in  the  data  once  the  MSS  is  sited  or 
operating.  At  that  time,  it  would  be  worthwhile  to  perform 
this  kind  of  analysis  to  improve  the  utility  of  the  data  from 
the  MSS. 

GENERALIZED  RAY  THEORY  ANALYSIS 

In  the  ray  tracing  analysis  presented  in  the  preceding  section, 
only  the  direct  arrivals  were  included.  No  reflected  or  diffracted 
phases  were  included  in  the  calculations.  This  represents  an 
appropriate  solution  since  the  first  arrivals  from  these  deeper 
events  will  be  the  direct  rays,  any  reflected  or  other  phases 
will  arrive  later  in  the  wavetrain.  However,  when  the  hypocentral 
depths  are  relatively  shallow  (less  than  ~100  km  for  this  velocity 
structure) ,  the  situation  changes.  The  direct,  diving  rays  are 
refracted  down  into  the  mantle  and  channeled  out  to  distances 
beyond  the  area  of  interest.  This  point  was  illustrated  in 
Figure  13.  Hence,  the  siting  area  sits  in  a  classical  shadow 
zone  for  the  diving  or  direct  rays  from  shallow  sources  in  the 
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Kuriles  or  Kamchatka.  The  first  arrivals  are  thus  multiple 
or  diffracted  Pg  and  Pn  type  phases.  In  order  to  examine  these 
phases,  it  has  been  necessary  to  use  a  generalized  ray  theory 
code  to  calculate  synthetic  seismograms  in  the  region. 

A  detailed  discussion  of  the  calculation  of  synthetic 
seismograms  is  beyond  the  scope  of  this  report.  The  methods 
involved  have  been  extensively  discussed  in  the  open  seismological 
literature  (e.g.,  de  Hoop,  1960;  Helmberger,  1968;  Helmberger 
and  Wiggins,  1971,  Fuchs,  1971;  Wiggins,  1976;  Chapman,  1976; 
Mellman  and  Helmberger,  1978).  In  our  analysis,  we  have  used 
the  asymptotic  Cagniard  formulation  of  Wiggins  and  Helmberger 
(1970) .  Use  of  this  method  provides  greater  insight  into  the 
physical  processes  involved  than  reflectivity  methods  (i.e., 

Fuchs,  1971). 

As  noted  earlier,  the  first  arrivals  from  shallow  sources, 

at  the  distances  of  interest  to  us  here,  are  dominated  by 

multiple  reflected  and  diffracted  P  and  P  type  phases.  The 

g  n 

resulting  seismogram  is  principally  dependent  upon  the  existence 
of  a  waveguide  layer  in  which  the  seismic  energy  can  be  efficiently 
trapped.  The  details  of  the  seismogram,  however,  depend  strongly 
upon  fine  structure  such  as  sediment  thicknesses,  on  the  depth 
of  the  source  re- ttive  to  the  waveguide  layer (s),  and  on  such 
factors  as  the  dt  h  of  burial  of  the  downhole  seismometer. 

It  was  not  ccupr.t -‘c.’ cnally  practical  or  cost  effective  to 
examine  many,  rou'h  less  all,  of  the  possible  source-receiver 
pairs  in  th*>  *uj ;  i.es/Kamchatka  area.  As  an  alternative  approach, 
we  have  used  first  an  illustrative  calculation  by  Helmberger 


and  Engen  (1980)  of  propagation  in  a  layer  over  a  half  space; 
and,  second,  we  have  calculated  the  response  of  our  upper  mantle 
model  for  several  different  shallow  source  depths.  In  each 
case,  profiles  of  the  synthetic  seismograms  are  presented  from 
100  km  to  1,000  km  for  each  of  the  three  fundamental  faults 
(vertical  strike-slip,  vertical  dip-slip,  45°  dip-slip) .  The 
first  example,  shown  in  Figures  32a,  32b,  and  32c,  is  the 
calculation  by  Helmberger  and  Engen  (1980) .  This  study  used 
the  following  models 


P  Velocity 

S  Velocity 

Density 

Thickness 

(km/sec) 

(km/sec) 

(gm/ee) 

(km) 

6.2 

3.5 

2.7 

32 

8.2 

4.5 

3.4 

00 

This  calculation  includes  all  of  the  multiples  and  trapped  energy 
within  the  waveguide  layer.  The  purpose  of  including  this 
result  here  is  to  illustrate  the  complexity  introduced  into 
the  initial  portions  of  the  seismogram  by  this  trapped  energy 
and  to  give  an  indication  of  the  expected  amplitude  fall  off 
with  distance. 

The  second  set  of  examples  utilized  in  our  study 
was  obtained  by  calculating  the  model  response  of  our  upper 
mantle  model  (Table  I)  for  three  source  depths;  5  km,  13  km, 
and  45  km.  This  calculation  differs  from  the  previous  example 
in  that  it  does  not  include  the  higher  multiples  in  the  upper 
waveguide.  These  calculations  include  all  of  the  direct  phases, 
head  waves,  pP  and  sP  phases  (off  of  both  the  ocean  bottom  and 
the  underside  of  the  Moho)  and  such  phases  as  s  P  (energy 
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Figures  32a,b,c.  Layer  over  half-space  calculation  by 

Helmberger  and  Engen  (1980) .  The  vertical  displacement 
responses  are  calculated  for  each  of  the  three  fundamental 
double  couple  systems:  strike-slip  (a),  vertical  dip-slip 
(b) ,  45°  dip  slip  (c) .  Each  response  has  been  convolved 
with  a  triangular  time  function  (0.5,  0.5  sec)  and  a  WWSSN 


short  period  response.  The  numbers  on  the  right  indicate 
the  trace  amplitudes  in  10“^cm.  Assuming  a  source  with 
static  moment  M  =4np  x  10^0  dyn-cm. 
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Figure  32b 
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TABLE  Z 

Velocity  Model  used  for  Generalized  Ray  Calculations 


Velocity 

(km/sec) 

S  Velocity 
(km/sec) 

Density 

(gm/ee) 

Thickness 

(km) 

1.5 

0.0 

1.0 

4 

6.2 

3.5 

2.7 

15 

8.2 

4.74 

3.4 

70 

7.5 

4.33 

3.16 

40 

7.7 

4.45 

3.23 

40 

7.8 

4.51 

3.26 

20 

8.0 

4.62 

3.32 

13.3 

8.05 

4.65 

3.34 

13.3 

8.11 

4.69 

3.36 

13.3 

8.16 

4.72 

3.38 

13.3 

8.22 

4.75 

3.40 

13.3 
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leaving  downward  from  the  source  as  shear  and  converting  to  P 
at  the  Moho) .  The  first  three  figures#  33,  34,  and  35#  show 
the  synthetic  seismograms  with  each  trace  amplitude  normalized 
independently.  The  last  figure,  36,  is  a  duplication  of 
Figure  33  with  the  trace  amplitudes  normalized  to  the  100  km 
record  in  order  to  more  dramatically  illustrate  the  shadow 
zone  effect  on  the  arrivals  in  the  seismograms .  In  each 
example,  the  first  arriving  phase  are  the  and  slightly  diving 
phases,  the  second  major  arrivals  are  the  Pg  direct  and  first 
multiple  bounce  phases. 
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Figures  33a, b,c.  Generalized  ray  calculation  for  model  shown  in 
Table  I  and  a  source  depth  of  5  km  beneath  the  ocean  bottom. 
Included  in  the  ray  sum  are  the  direct  ray,  down  going  primary  rays 
on®  reflection,  no  mode  conversion),  and  primary  rays  from 
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Figure  34c 


SOURCE  R1  45.  KM 

_  14.  116 


STRIKE 

100.  KM 
T-  13.  8 


200.  KM 
T*  25.5 


SLIP 

il - It 


3.  015 


300.  KM 
T-  37.  5' 


i .  3  1 5 


400.  KM 
T-  49.  S 


500.  KM 
T*  61.7 


600.  KM 
T-  73.  8 


700.  KM 
T  *  85.  9 


0.  966 


G.  637 


0.  366 


0.  453 


800.  KM 
f  =  98.  C 


l/\/^A, - — 


Q.  475 


900.  KM 
T  =  1  10.  1 


tj| 

J  j  i/'\/wV' — r'J\p- 


0.  458 


0.  398 


Figures  35a, b,c.  Generalized  ray  calculations  assuming  a 
source  depth  of  45  km.  Also  included  in  the  ray  sum  are 
the  first  reflections  (sP,  pP)  from  the  Moho. 
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Figures  36a, b,c.  This  figure  is  a  u  plot  of  Figure  33 
amplitudes  scaled  the  same.  The  decay  in  amplitude  is 
zone"  caused  by  the  low  velocity  zone  at  85  km  depth. 
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SUMMARY  AND  CONCLUSIONS 

This  project  has  involved  the  calculation  of  literally  tens 
of  thousands  of  rays  through  a  rapidly  changing,  laterally  heter- 
geneous  material.  The  purpose  of  the  study  was  to  provide  an 
additional  seismological  evaluation  of  potential  sites  for  the 
MSS  system.  In  the  course  of  the  investigation,  a  detailed 
three-dimensional  upper  mantle  model  was  constructed  for  the 
Kuriles/Kamchatka  area.  The  adoption  of  this,  or  any  other, 
model  is  a  crucial  step  in  the  overall  analysis.  The  results 
presented  within  this  report  are  based  very  strongly  upon  the 
velocity  structure  selected.  The  model  used  here  was  a  mod¬ 
ification  of  the  Veith  (1974)  model.  This  model  has  the  ad¬ 
vantage  of  being  based  on  a  very  large  set  of  corrected  seismic 
travel-times  and  is  felt  to  be  quite  reliable.  Large  changes  in 
the  velocity  structure  used  could  alter  some  of  our  conclusions. 
However,  we  have  also  performed  the  simulation  using  an  unwarped 
slab  model.  The  results  do  not  greatly  differ  between  the  two 
calculations.  Thus,  we  can  analyze  the  results  with  some 
confidence. 

There  is  a  tendency  to  become  somewhat  inundated  with  detail 
when  one  attempts  to  assimilate  the  results  from  this  study.  It 
is  necessary  to  maintain  a  broad  overview  and  to  identify  con¬ 
sistent  patterns  in  the  predicted  energy  density.  This  has  been 
true  throughout  the  project.  The  scale  of  the  questions  addressed 
and  the  broad  geographical  area  involved  required  that  the  in¬ 
vestigation  emphasize  broadly-based  analysis  methods  and  not 
concentrate  on  small  details.  On  this  basis,  we  can  draw  several 


74 


significant  conclusions  from  the  study  and  provide  a  firm 
foundation  for  siting  recommendations. 

In  referring  to  the  four  potential  MSS  sites  identified 
here,  we  shall  use  the  designators  Sta.  3,  Sta.  4,  Sta.  1, 
and  Sta.  la  to  indicate  the  four  sites  in  order  of  decreasing 
latitude  (i.e.,  Sta.  3  is  the  northernmost  site,  Sta.  la  the 
southernmost) . 

When  6ne  considers  the  results  of  the  ray-tracing  analysis 
(sources  deeper  than  100  km),  we  find  that  Sta.  3  is,  by  a  large 
margin,  the  best  site  in  terms  of  ability  to  see  diving  rays 
from  most  of  the  region.  Only  three  of  the  sixteen  hypocenters 
considered  resulted  in  Sta.  3  being  in  a  shadow  zone  (these 
hypocenters  were  all  at  150  km  depth  for  the  three  northernmost 
source  sets).  Using  the  same  criteria,  Sta.  1  represents  the 
least  desirable  site.  However,  if  we  are  particularly  interested 
in  shallow  events  under  Kamchatka  the  two  southern  sites,  Sta.  1 
and  especially  Sta.  la  are  superior  to  the  other  stations  in 
their  ability  to  see  diving  rays  from  northern  events.  Indeed, 
Sta.  la  is  at  a  very  intense,  bright  spot  for  several  northern 
locations.  In  the  end,  the  seismological  evaluation  of  each  site 
is  dependent  upon  a  decision  of  which  seismic  area  is  of  greatest 
interest.  From  our  perspective,  for  best  uniform  coverage  of 
the  region,  we  would  rate  the  four  site  locations  in  order  of 
desirability  as:  Sta.  3,  Sta.  la,  Sta.  4,  and  finally,  Sta.  1. 

On  the  basis  of  other  siting  criteria,  a  station  that  is 
situated  in  a  shadow  zone  for  much  of  the  studied  seismic  area 
may  be  desirable.  The  utility  of  such  a  station  would  necessarily 


require  a  thorough  understanding  of  the  multiple  phases  and  non¬ 
geometric  arrivals  that  typify  the  location.  A  simulation  and 
discussion  of  these  multiple  reflected  and  refracted  phases  was 
presented  earlier  in  this  report,  Figures  32-36.  Distances  from 
the  source  as  well  as  fine  details  of  the  path  are  of  utmost 
importance.  The  obvious  difficulty  with  the  interpretation  of 
such  seismic  data  for  the  determination  of  source  parameters 
is  the  poorly  defined  path-dependent  propagational  effects. 

There  are  several  obvious  directions  for  future  research 
that  have  come  to  light  during  the  course  of  this  investigation. 
Ultimately  synthetic  seismograms  that  correctly  account  for 
laterally  varying  structures  will  be  required.  With  the  addition 
of  multiple  reflected  arrivals,  the  current  code  could  be 
modified  to  generage  synthetic  seismograms  that  correctly 
account  for  amplitude,  timing  and  phase-shifts  of  later  arrivals. 
Siu  h  a. code  would  be  quite  valuable  for  the  determination  of 
source  characteristics. 
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